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Abstract. Species composition of Eurasian water frogs, genus Pelophylax, in Lake Skadar (Montenegro) was analysed using 
mitochondrial (mtDNA) and nuclear (nuDNA) markers. Specimens were characterised at first using mtDNA sequences 
of NADH dehydrogenase subunit 3 (ND3) gene. Based on their mitochondrial genomes, 49 specimens were determined 
as Pelophylax kurtmuelleri while 39 specimens were identified as Pelophylax shqipericus. The systematic affiliation was 
evaluated further using serum albumin intron-i (SAI-i) nuDNA. The results of SAI-i analyses confirmed identification of 
both species but also their hybrids. The SAI-r variant of P. shqipericus was ~6oo nucleotide base pairs longer compared to 
P. kurtmuelleri. Five specimens contained both variants, indicating their hybrid origin. However, population allotment of 
hybrids was low, suggesting normal Mendelian inheritance in the interspecific mating of P. kurtmuelleri and P. shqipericus 
rather than a hybridogenetic mode of reproduction. Mito-nuclear discordance as a result of backcross hybridization was 
not observed in this study. However, the absence of mito-nuclear discordance is not surprising since backcross specimens 
were rarely detected in previous research of water frogs from Lake Skadar. This study showed that species composition in 
Lake Skadar is the same compared to 30 years ago, when a previous study analysed the water frogs using protein electro¬ 
phoresis. 
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Introduction 

With a surface area of -44,000 km 2 , the Southeast Adri¬ 
atic Drainages sensu Abell et al. (2008) (hereinafter re¬ 
ferred to as SEAD) is one of the smallest freshwater ecore- 
gions in Europe (Fig. 1). The SEAD was delineated from 
its neighbouring freshwater ecoregions by endemic-rich 
fauna of fish (Abell et al. 2008, Skoulikidis et al. 2009), 
gastropods (Albrecht & Wilke 2009), crustaceans (Wy- 
socka et al. 2014) and amphibians with particular refer¬ 
ence to water frogs (Dzukic & Kalezic 2004). The ob¬ 
served independent evolutionary histories of distinct biota 
in the SEAD were discussed in the light of complex geo- 
morphologic events, climate history, and ecosystem diver¬ 
sity (Stankovic i960, Albrecht & Wilke 2009). 

In this study, we analysed species composition of wa¬ 
ter frogs in Lake Skadar, the largest freshwater body in 
the SEAD (surface area of -370 km 2 ), as well in all of 
entire Southeastern Europe. Lake Skadar was formed 


due to complex tectonic folding and faulting within the 
northeast wing of Old Montenegro anticlinorium (High 
Karst Zone) during the Neogene or even Paleogene pe¬ 
riod (Pesic & Gloer 2013). Together with Lake Ohrid 
and Lake Prespa, Lake Skadar is considered a hot-spot 
of endemic flora and fauna in Southeastern Europe (Al¬ 
brecht & Wilke 2009, Pesic et al. 2009). Two Eura¬ 
sian water frog species have been described from SEAD, 
the Albanian pool frog, Pelophylax shqipericus (Hotz, 
Uzzell, Gunther, Tunner & Heppich, 1987), and the 
Balkan pool frog, Pelophylax kurtmuelleri (Gayda, 1940) 
(Schneider et al. 1993, Uzzell et al. 2009, AmphibiaWeb 
2016). Both species are common in Lake Skadar (Uzzell 
et al. 2009) and their type localities are shown in Figure 1. 
The third species reported from the SEAD is the Epirus 
water frog, Pelophylax epeiroticus (Schneider, Sofiani- 
dou & Kyriakopoulou-Sklavounou, 1984). But it has 
never been observed in Lake Skadar (Uzzell et al. 2009, 
AmphibiaWeb 2016). 
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Pelophylax shqipericus was first biochemically recog¬ 
nised and described as an unnamed species from Lake Ska¬ 
dar by Hotz & Uzzell (1982). Later Schneider & Haxhiu 
(1994) had difficulties in distinguishing the mating calls of 
P. shqipericus from that of the pool frog P. lessonae (Came- 
rano, 1882). However, full-species status of P. shqipericus 
was supported by numerous studies on morphology, allo- 
zymes (Hotz & Uzzell 1982, Beerli et al. 1996), mtDNA 
(Plotner 1998, Akin et al. 2010a, 2010b, Plotner et al. 
2010) and nuDNA sequences (Pyron & Wiens 2011). Al¬ 
though the species status of P. shqipericus was confirmed, 
its hierarchical relationships among the Eurasian wa¬ 
ter frogs remained less clear; P. shqipericus was placed in 
a monophyletic group with P. lessonae and P. bergeri in 
mtDNA phylogenies (e.g. Plotner et al. 2010), while its 
phylogenetic position based on nuDNA markers was un¬ 
certain (e.g. Beerli et al. 1996, Pyron & Wiens 2011). 

In Lake Skadar, P. shqipericus lives in sympatry with 
P. kurtmuelleri (Schneider et al. 1993, Uzzell et al. 2009). 
In comparison with P. shqipericus, the systematic status of 
P. kurtmuelleri is less clear and its full-species status has of¬ 
ten been debated. Crochet & Dubois (2004) and Spey- 
broeck et al. (2010) argued that there is low genetic diver¬ 
gence between P. kurtmuelleri and P. ridibundus, implying 
that P. kurtmuelleri is conspecific with P. ridibundus. Also, 
the results of some allozyme analyses were controversial in 



type locality 


Figure 1. Map presenting sampling sites of water frogs in Lake 
Skadar and type localities of species occurring in Southeast 
Adriatic drainages (SEAD) sensu Abell et al. ( 2008 ). National 
boundaries are indicated by black lines while the border of SEAD 
is shown by a thick orange line. 


distinguishing P. kurtmuelleri from P. ridibundus (Beerli et 
al. 1996, Hotz et al. 2013). In earlier studies of water frogs 
from Lake Skadar, specimens distinct from P. shqipericus 
were assigned to the Adriatic populations of the Eurasian 
marsh frog P. ridibundus (Pallas, 1771) (Hotz & Uzzell 
1982, 1983, Hotz et al. 1985, Guerrini et al. 1997, Spasic- 
Boskovic et al. 1999). But within P. ridibundus in Greece, 
two taxa were discriminated based on their mating calls 
and morphometry (Schneider & Sinsch 1992, Schnei¬ 
der et al. 1993). The new taxon was named Rana balcanica 
Schneider & Sinsch, 1992. Later this taxon was recognised 
as a junior synonym of P. kurtmuelleri (Dubois & Ohler 
1994). Further evidence for full-species status of P. kurtmu¬ 
elleri was found in its non-hybridogenetic interactions with 
P. ridibundus and P. lessonae (Hotz et al. 1985, Berger et al. 
1994), morphometries (Gavrilovic et al. 1999), bioacous¬ 
tics (Schneider & Sinsch 1992, Schneider et al. 1993, Lu- 
kanov et al. 2015), electrophoretic investigation (Soflani- 
dou et al. 1994), phylogenetic inference based on DNA se¬ 
quence data (Lymberakis et al. 2007, Pyron & Wiens 2011, 
Plotner et al. 2012, Hofman et al. 2015), and the results of 
research on cytogenetic differences and centromeric hybrid¬ 
ization (Marracci et al. 2011). Pelophylax kurtmuelleri is 
distributed up to 1,000 m above sea level throughout much 
of Greece and Albania (Uzzell et al. 2009). Furthermore, it 
is considered a naturalised alien species in Italy from where 
it spread to Slovenia (Bressi 2006). Also, P. kurtmuelleri is 
an introduced species in Denmark (Lever 2003). 

In areas where two species of water frog are in con¬ 
tact, their hybrids are extant (Arnold & Ovenden 2002, 
Mayer et al. 2013). Among hybrid forms in Eurasian wa¬ 
ter frogs, there are three different fertile hybrids with re¬ 
productive mode termed hybridogenesis: the Edible frog 
Pelophylax kl. esculentus (Linnaeus, 1758), the Graf’s hy¬ 
brid frog Pelophylax kl. graft (Crochet, Dubois, Ohler 
& Tunner, 1995) and the Italian edible frog Pelophylax kl. 
hispanicus (Bonaparte, 1839) (reviewed in Graf & Polls 
Pelaz 1989). Diploid hybridogenetic or hemiclonal hy¬ 
brids exclude one of the parental genomes from the hybrid 
germline, pass the remaining one clonally to their gam¬ 
etes, and provide the excluded genome via backcross mat¬ 
ing with the parental species (Schultz 1969, Schmeller 
et al. 2001). The genome of P. ridibundus is considered to 
induce hybridogenetic reproduction in all three hybrido¬ 
genetic complexes in Europe (Graf & Polls Pelaz 1989, 
Holsbeek & Jooris 2010). To date, the species of Eurasian 
water frogs and their hybrid forms were intensively stud¬ 
ied as models for understanding speciation in the context 
of hybridization and polyploidization (Pagano et al. 2001, 
Plenet et al. 2005, Christiansen & Reyer 2009, Chris¬ 
tiansen et al. 2010, Jakob et al. 2010, Hoffman & Reyer 
2013, Pruvost et al. 2013, Hoffman et al. 2015). 

During the 1980s, water frogs from Lake Skadar were re¬ 
peatedly studied using protein (allozyme) electrophoresis 
and chromosome analyses to investigate interspecific hy¬ 
bridization (Hotz & Uzzell 1982,1983, Hotz et al. 1985, 
Guerrini et al. 1997). The presence of P. shqipericus (at 
that time referred to as an unknown species), P. kurtmuel- 
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leri (referred to as P. ridibundus from the Adriatic region) 
and their hybrids was reported for Lake Skadar (Hotz & 
Uzzell 1982,1983). The comparison of results obtained by 
enzyme electrophoresis from somatic tissues with those 
from individual primary oocytes indicated that P. shqiperi¬ 
cus does not produce hybridogenetic hybrids in interspe¬ 
cific crosses with P. kurtmuelleri (referred to as “Adriatic 
P. ridibundus ”) (Hotz et al. 1985). The same was observed 
for hybrids of P. kurtmuelleri (referred to as “Balkan Rana 
ridibunda”) and P. epeiroticus (Rana epeirotica) (Guerrini 
et al. 1997). Reported frequency of hybrids of P. kurtmuel¬ 
leri and P. shqipericus from Lake Skadar was low (8.7%), 
unlike the abundant hemiclonal hybrid lineage P. kl. escu- 
lentus in Central and Western Europe (Hotz & Uzzell 
1982). It was indicated that the lack of hybridogenetic re¬ 
production in these hybrids was caused by a failure of hy¬ 
brids to “induce” the exclusion of one parental species ge¬ 
nome during gametogenesis of hybrids, or a “resistance” of 
P. shqipericus genome to such exclusion, or both (Guer¬ 
rini et al. 1997). Nevertheless, it was acknowledged that the 
P. kurtmuelleri genome in hybrids from Lake Skadar does 
not contain putative inducing factors for the exclusion of 
P. shqipericus genome during gametogenesis (Guerrini 
et al. 1997). However, a few backcross specimens were ob¬ 
served in Lake Skadar, indicating that non-hybridogenetic 
Fi hybrids of both parental combinations were not com¬ 
pletely sterile (Hotz & Uzzell 1982,1983). 

In recent years, analyses of DNA sequence fragments 
have become common in water frog research. Mitochon¬ 
drial DNA sequences such as cytochrome b, 16S rRNA 
(Lymberakis et al. 2007), NADH dehydrogenase subunit 
2 gene (ND2) and NADH dehydrogenase subunit 3 gene 
(ND3) (Plotner et al. 2008, Akin et al. 2010a, b, Do- 
meneghetti et al. 2013, Hotz et al. 2013) were proven in¬ 
formative in phylogenetic inference, as well as nuDNA se¬ 
quences, e.g. serum albumin intron-i (SAI-i) (Plotner 
et al. 2009, Hauswaldt et al. 2012). Also, the difference 
in lengths of the SAI-i region has become a widely used 
method for distinguishing among species and hybrids 
within Pelophylax species (Hauswaldt et al. 2012, Mayer 
et al. 2013, Herczeg et al. 2016). Therefore, one of the main 
focuses of this research was to obtain sequences of the 
SAI-i nuclear region of P. shqipericus from its type local¬ 
ity. Further, our aim was also to apply nuclear (SAI-i) and 
mtDNA (ND3) markers to investigate, whether the species 
composition in Lake Skadar had changed compared to the 
study conducted 30 years ago. Also, we investigated occur¬ 
rence of mtDNA introgression between P. shqipericus and 
P. kurtmuelleri via backcross mating of their hybrids. 

Materials and methods 

Sampling procedure 

Eighty-eight water frogs were collected at four localities 
on Lake Skadar in Montenegro (Supplementary Table Si, 
Fig. 1): 27 specimens from Rijeka Crnojevica (42.35^, 
19.04°E); 29 from Mareza (42.45 0 N, 19.20°E); 26 from 


Lake Skadar (42.33^, 19.09^); 6 from Moraca (42.28°N, 
19.13°E). The tip of a toe was clipped from each frog and 
preserved in 96% ethanol for further processing. Frogs were 
released after this procedure, no vouchers were maintained. 

DNA extraction and PCR amplifications 

Total genomic DNA was extracted from toe tissue using 
GenElute Mammalian Genomic DNA Miniprep Kit (Sig- 
ma-Aldrich) following the manufacturer’s instructions. 
DNA quality and concentration was checked with a Na¬ 
nodrop ND-100 spectrophotometer (Nanodrop Technol¬ 
ogies). The primer pair used for ND3 was ND3L/ND3H 
(Akin et al. 2010a) and the primer pair for SAI-i was Pel- 
SA-Fi/Pel-SA-R2 (Hauswaldt et al. 2012). Each PCR was 
conducted using HotStarTaq Plus Master Mix (Qiagen) in 
a total volume of 10 pL containing 0.5 U of HotStarTaq Plus 
DNA Polymerase, 200 pM of each dNTP, 1.5 mM MgCl 2 , 
0.25 pM of each primer, and 40-120 ng of template DNA. 
Thermocycling profiles started with a 5 min activation 
step at 95°C, and were followed by 35 cycles of denatura- 
tion (30 s at 94°C), annealing (20 s at 62°C for ND3 or 40 s 
at 59°C for SAI-i) and extension (60 s at 72°C for ND3 or 
100 s at 72°C for SAI-i), and a final extension step at 72°C 
for 10 min. PCR amplicons were purified using MinElute 
PCR Purification Kit (Qiagen). 

Visualisation of PCR amplicons using gel electrophoresis, 
sequencing and cloning procedures 

Three microliters of PCR amplicons were run on a 1.5% aga¬ 
rose gel. Eighty-eight purified amplicons of ND3 and six¬ 
teen randomly chosen SAI-i amplicons (eight per species) 
were sequenced in both directions using the same primers 
used in PCR (Supplementary Table Si). Sequencing reac¬ 
tions were conducted on ABI 3730XL DNA Analyzer (Ap¬ 
plied Biosystems) in Macrogen Europe sequencing service 
(Amsterdam/NL). Chromatographs were checked and ed¬ 
ited using SEQUENCHER version 5.3 (Gene Codes Corp.). 
SAI-i amplicons were cloned in case of heterozygous in¬ 
dividuals. The PCR amplicons of heterozygous individu¬ 
als were ligated into vectors and transformed into bacteria 
using the pGEM-T Vector System II (Promega). Plasmid 
DNA from nine positive clones per specimen was isolated 
and purified using the Wizard Plus SV Miniprep DNA Pu¬ 
rification System (Promega), and inserts were sequenced 
with the primers used during PCR and new internal prim¬ 
ers (INTSHQ-Fi: 5’-GCACGGAGCATCGATAGTTT-3’, 
and INTSHQ-Ri: 5’-AGGCATAAGGTGCCCATCTA-3’). 


Phylogenetic inference 
using ND 3 and SAI -1 sequences 

Sequence traces were analysed using SEQUENCHER ver¬ 
sion 5.3 (Gene Codes Corp.). The SAI-i sequences gen- 
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Table 1. List of unique sequences for nuDNA SAI-1 region in water frogs from Lake Skadar. 


Sample 

ID 

Locality 

ID of SA-1 
haplotypes 

Taxonomic 

affiliation 

NCBI GenBank 
Acc. No. of SAI-1 

RC2 

Rijeka Crnojevica 

Seql 

P. kurtmuelleri 

MH038010 

RC3 

Rijeka Crnojevica 

Seql 

P. kurtmuelleri 

MH038011 

RC7 

Rijeka Crnojevica 

Seql 

P. kurtmuelleri 

MH038012 

RC18 

Rijeka Crnojevica 

Seql 

P. kurtmuelleri 

MH038013 

MA19 

Mareza 

Seql 

P. kurtmuelleri 

MH038016 

MA21 

Mareza 

Seql 

P. kurtmuelleri 

MH038017 

MA22 

Mareza 

Seql 

P. kurtmuelleri 

MH038018 

MA26 

Mareza 

Seql 

P. kurtmuelleri 

MH038019 

RC21 

Rijeka Crnojevica 

Seq2, Seq3 

P. shqipericus 

MH038014 

RC23 

Rijeka Crnojevica 

Seq2, Seq4 

P. shqipericus 

MH038015 


erated in this study were examined using Nucleotide Ba¬ 
sic Local Alignment Search Tool (Nucleotide BLAST: 
http://blast.ncbi.nlm.nih.gov/Blast.cgi) to screen for the 
most similar sequences in the GenBank nucleotide da¬ 
tabase (National Center for Biotechnology Information, 
U.S. National Library of Medicine, USA). Sequence align¬ 
ments were inferred using MAFFT version 7.187 (Katoh 
et al. 2005) on the CIPRES Science Gateway version 3.1 
of Miller et al. (2010) (http://www.phylo.org). Phyloge¬ 
netic tree reconstructions were conducted using sequenc¬ 
es obtained in this study and the sequences of ND3 and 
SAI-i which were downloaded from the GenBank (https:// 
www.ncbi.nlm.nih.gov/genbank/) and are those of Sumida 
et al. (2001), Plotner et al. (2008, 2009, 2010), Akin et 
al. (2010a, b), Hofman et al. (2012, 2015), Domeneghet- 
ti et al. (2013), Hotz et al. (2013), Dubey et al. (2014), 
MikulIcek et al. (2014); see Supplementary Table Si. The 
sequences obtained in this study were likewise deposited 
in GenBank; see Supplementary Table Si and Table 1. Phy¬ 
logenetic inference was performed using Maximum Like¬ 
lihood (ML) and Bayesian Analysis (BA). The best-fit evo¬ 
lutionary model used in ML and BA was calculated using 
the Bayesian Information Criterion (BIC) in jModelTest2 
version 2.1.6 (Darriba et al. 2012) as implemented on the 
CIPRES. Best-fit model of nucleotide substitution for ND3 
gene was TrN (Tamura & Nei 1993) while the model for 
SAI-i sequences was HKY (Hasegawa et al. 1985) with a 
gamma-distributed rate variation among sites and a signif¬ 
icant proportion of invariable sites. 

The ML analyses for ND3 and SAI-i were run using 
MEGA version 6 (Tamura et al. 2013) with optimised pa¬ 
rameters. Initial tree(s) for the heuristic search were ob¬ 
tained by applying the Neighbor-Joining (NJ) method to a 
matrix of pairwise distances estimated using the Maximum 
Composite Likelihood approach. The obtained phylogeny 
was tested using 1,000 nonparametric bootstrap replicates. 

In the analyses of ND3, the BA was run in MrBayes 3.2 
(Ronquist et al. 2012) on CIPRES. Two independent runs 
with four MCMC chains were run for 50 million genera¬ 
tions and sampled every 5,000 generations, with temper¬ 
ature parameter set to 0.2 and the first 12.5 million gen¬ 


erations discarded as burn-in. In the analyses of SAI-i se¬ 
quences, the BA was conducted using BAly-Phy version 
2.3.8 to incorporate insertion/deletion information into 
phylogeny estimation (Redelings & Suchard 2005, 2007, 
Suchard & Redelings 2006). BA in BAly-Phy was run 
using 100,000 iterations, the HKY model of nucleotide 
substitution, and RS07 insertion/deletion model (Rede¬ 
lings & Suchard 2007). To reduce computational time in 
BAly-Phy, SAI-i sequence data set was reduced to one ran¬ 
domly chosen haplotype per species. 

The convergence of runs in BA was screened using 
AWTY (Nylander et al. 2008) while effective sample siz¬ 
es of parameters were checked using TRACER 1.5 (Drum¬ 
mond & Rambaut 2007). Nodes in phylograms with boot¬ 
strap values P > 70 in ML and posterior probabilities (pp) 
> 0.95 in BA were considered as support. 


Results 

Phylogenetic inference based on mtDNA (ND 3 ) 

ND3 sequences were generated for all individuals (Sup¬ 
plementary Table Si). Forty-nine sequences which corre¬ 
sponded to seven unique haplotypes were placed in the 
P. kurtmuelleri lineage. Thirty-nine sequences with eleven 
haplotypes were positioned in the P. shqipericus lineage in 
ML and BA phylograms (Fig. 2). Our phylogenetic infer¬ 
ence based on ND3 placed P. kurtmuelleri haplotypes in 
the group “ridibundus-bedriagae” sensu Plotner & Ohst 
(2001), but they were genetically distinct from P. ridibun- 
dus (Fig. 2). Haplotypes of P. shqipericus were positioned 
within the major group “lessonae” sensu Plotner & Ohst 
(2001), and formed a distinct supported group in both ML 
and BA phylogenies using ND3 mtDNA. 

Length variation of SAI -1 region in water frogs 
from Lake Skadar 

The PCR amplification of SAI-i fragments and subse¬ 
quent visualisation on agarose gels were completed with 
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70/- 


80/0.99 

70/0.92 


93/1 r C 

r ic 

L|-Cil 


99/1 


85/1- 

96/1 

-/0.98 


P. esculentus (KJ160671: Slovakia) 

P. ridibundus( GU812152: Russia) 

P. rid/bundus (HE855847: Bulgaria) 

P. ridibundus (HE855823: Turkey) 

P. ridibundus (HE855774, HE855813, HE855817: Greece) 
P. ridibundus (KJ152758: Switzerland) 

- P. ridibundus (HE855751, HE855769, GU812151: Greece) 
P. ridibundus/P. escu/enfus (Haplotypel: *) 

— P. kurtmuelleri (G\J8'\2'\ 57: Turkey) 

P. kurtmuelleri (GU812156: Russia) 
j P. kurtmue//e/7(Haplotype6: Greece) 

L p. kurtmuelleri (KM 52759: Switzerland) 

- ♦ P. kurtmue//er/(Skadarl 7: Montenegro) 

- P. kurtmuelleri (AM749702, AM749703: Greece) 

- ♦ P. kurtmuelleri/hybrid (Skadar7: Montenegro) 

- ♦ P. kurtmuelleri (Skadarl: Montenegro) 

- P. kurtmuelleri (AM749697: Greece) 

P. kurtmuelleri (NC026895: Greece) 

- ♦ P. kurtmuelleri (Skadarl 6: Montenegro) 

♦ P. kurtmuelleri / hybrid (Skadar2: Montenegro) 

♦ P. kurtmuelleri (Skadarl 8: Montenegro) 

P. kurtmuelleri (Haplotype7: Greece, Egypt) 

P. kurtmuelleri (HE855780: Greece) 

P. kurtmuelleri (GU812155: Russia) 

P. kurtmuelleri (HE855772, GU812154: Greece) 

P. kurtmuelleri (HE855748: Greece) 

P. kurtmuelleri (HE855744, GU812153: Greece) 

P. kurtmuelleri (H a ploty pe5: ****) 

P. kurtmuelleri (KJ 1 ] 52760: Switzerland) 

♦ P. kurtmuelleri (Skadar9: Montenegro) 

L P. kurtmuelleri (HE855739: Greece) 

[r P- kurtmuelleri (H E855741, AM749705: Greece) 

■ P. kurtmuelleri (AM749704, AM749706: Albania, Greece) 

99/1 ^The Levant MHG2 sensu Akin et al. 2010 (GU812158: Jordan) 


The Levant MHG2 sensu Akin et al. 2010 (GU812159: Jordan) 
Cilician west MHG4 sensu Akin etal. 2010 (GU812164: Turkey) 
Cilician west MHG4 sensu Akin et al. 2010 (GU812167: Turkey) 
Cilician east MHG5 sensu Akin et al. 2010 (GU812168: Turkey) 
Cilician east MHG5 sensu Akin etal. 2010 (GU812170: Turkey) 
Middle east MHG9 sensu Akin et al. 2010 (GU812230: Iran) 

— Central Asia MHG8 sensu Akin etal. 2010 (GU812229: Kazakhstan) 
gg^ | P. cypriensis (NC026893: Greece) 


78/- 


L Cyprus MHG3 sensu Akin et al. 2010 (GU812160: Cyprus) 

■ Cyprus MHG3 sensu Akin et al. 2010 (GU812161: Cyprus) 

95/1Central Asia MHG7 sensu Akin etal. 2010 (GU812225: Iran) 

■— Central Asia MHG7 sensu Akin etal. 2010 (GU812228: Iran) 

K ef. caralitanus MHG6a sensu Akin etal. 2010 (GQ902080: Turkey) 
cf. caralitanus MHG6a sensu Akin etal. 2010 (GQ902085: Turkey) 

___ _j cf. cerigensis MHG6b sensu Akin et al. 2010 (GU812181: Greece) 

P cf. cerigensis MHG6b sensu Akin etal. 2010 (GU812182: Turkey) 

- cf. bedriagae s.s. MHG6c sensu Akin et al. 2010 (GU812189: Turkey) 
cf. bedriagae s.s. MHG6c sensu Akin et al. 2010 (GU812192: Turkey) 

... i— Euphrates MHG6d sensu Akin et al. 2010 (GU812212: Turkey) 

»- Euphrates MHG6d sensu Akin et al. 2010 (GU812217: Iran) 

P. epeiroticus (NC026894: Albania) 


95/1 


- P. cretensis (NC025575: Greece) 


96/1 \ 


75/0.95.. 

96/1..;: 

84/- 


0.05 


P. lessonae (KJ152766: Switzerland) 

- P. lessonae (AM749716: Estonia) 

P. lessonae (KJ152765: Switzerland) 

P. lessonae/P. ridibundus!P. esculentus (Haplotype3: ***) 

P. lessonae/P. ridibundus/P. esculentus (Haplotype2: **) 

P. lessonae (KJ 152767: Switzerland) 

L P. ridibundus/P. esculentus (Haplotype8: Sweden) 

P. lessonae (KJ152764: Switzerland) 

P. bergeri (GU812231: Italy) 

P. bergen(KJ152763: Switzerland) 

P. bergen(Haplotype4: Italy, France) 

— ♦ P. shqipericus (Skadarl 0: Montenegro) 

♦ hybrid (Skadarl 5: Montenegro) 
j ♦ P. shqipericus (Skadar5: Montenegro) 

L ♦ P. shqipericus (Skadarl 3: Montenegro) 

♦ P. shqipericus (Skadar8: Montenegro) 

P. shqipericus (NC026896: Greece) 

♦ P. shqipericus (Skadarl 2: Montenegro) 

♦ P. shqipericus (Skadarl 4: Montenegro) 

♦ P. shqipericus (Skadar6: Montenegro) 

♦ P. shqipericus (Skadar3: Montenegro) 

♦ P. shqipericus (Skadar4: Montenegro) 

♦ P. shqipericus (Skadarl 1: Montenegro) 
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Figure 2. Phylogenetic tree inferred by ML analysis using ND3 sequences of water frogs. Newly obtained haplotypes in this study 
were marked by black rhombi. Node supports are given as bootstrap values (P) in ML analysis (showing values > 70) and posterior 
probabilities (pp) in BA (showing values > 0.9). Legend: * - Bulgaria, France, Greece, Flungary, Italy, Latvia, Poland, Romania, Russia, 
Serbia and Slovakia; ** - Germany, Poland, Romania and Slovakia; *** - Denmark, Estonia, Germany, Italy, Latvia, Lithuania, Poland, 
Romania, Slovakia and Sweden; **** - France (Acc. No. AM749707), Germany (Acc. No. AM900647), Latvia (Acc. No. AM900652), 
Lithuania (Acc. No. AM900648), Romania (Acc. No. AM900650), Ukraine (Acc. No. AM900651), Italy, Greece and Macedonia. 
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all samples (Supplementary Table Si, Fig. 3). Two am- 
plicon lengths (-720 bp and -1200 bp) were detected. 
Forty-six specimens were determined as P. kurtmuelleri, 
since they carried SAI-i variants of -720 bp length and 
P. kurtmuelleri mtDNA genomes. Thirty-seven specimens 
were identified as P. shqipericus, since they carried SAI-i 
variants of -1200 bp and P. shqipericus mtDNA genomes. 
Furthermore, a mixture of both SAI-i variants was ob¬ 
served in five samples, indicating hybridization between 
P. kurtmuelleri and P. shqipericus in Lake Skadar. Three 
hybrid specimens had P. kurtmuelleri mtDNA genomes 
(haplotypes Skadarn and Skadar) while two hybrids had 
P. shqipericus mtDNA genomes (Skadaris). There were 
no cases of mito-nuclear discordance (mtDNA introgres- 
sion) in P. kurtmuelleri or P. shqipericus from Lake Ska¬ 
dar. 


Phylogenetic inference using nuDNA SAI -1 region 

Sequencing of the SAI-i region was done on 16 random¬ 
ly chosen samples originating from P. kurtmuelleri and 
P. shqipericus while hybrid specimens were avoided. We 
managed to get ‘good’ SAI-i sequences from eight speci¬ 
mens of P. kurtmuelleri and all sequences were identical 
(variant Seqi) (Table 1). This haplotype was identical to the 
GenBank sequence which originates from P. kurtmuelleri 
from Greece (FN432367). 

However, we failed to obtain ‘good’ sequences for the 
eight samples of P. shqipericus. Therefore, we cloned am- 
plicons of two individuals and sequenced nine clones per 
amplicon. Among the 18 sequences, we obtained three 
variants of SAI-i (Table 1). After BLAST search, the vari¬ 
ant Seq2 showed 99% similarity with the SAI-i sequence 
of P. shqipericus obtained by Dubey & Dufresnes (2017). 

Additionally, we reconstructed phylogenetic trees us¬ 
ing the SAI-i sequences of water frogs and new P. shqiperi¬ 
cus sequences from this study (Figs. 4 and 5). The SAI-i 
sequences of P. shqipericus form a monophyletic group 


(Fig. 4). Hierarchical relationships oiP shqipericus towards 
other species are not resolved in ML and BA phylogenies 
(Figs. 4 and 5). However, a long branch of P. shqipericus 
group points to a deep evolutionary divergence between 
this species and its closest relatives. 

Discussion 

Phylogenetic inference of P. shqipericus 
and P. kurtmuelleri from Lake Skadar 

In this study, we confirmed the presence of P. shqipericus 
and P. kurtmuelleri in Lake Skadar, using mitochondrial 
and nuclear genetic markers. Occurrence of both these taxa 
in Lake Skadar was first recognised by Hotz et al. (1985). 
Our results suggest that species composition in Lake Ska¬ 
dar has remained unchanged since analysed by Hotz et al. 
(1985) using protein electrophoresis. 

The ND3 haplotypes of P. shqipericus were positioned 
within the group “lessonae” sensu Plotner & Ohst 
(2001), comprising three species, P. lessonae, P. bergeri and 
P. shqipericus (Fig. 1). However this group is designated as 
monophyletic only in the ML analysis and is not support¬ 
ed as monophyletic in the SAI-i phylogenies (Figs. 4 and 
5). Nevertheless, a tendency of P. lessonae and P. shqiperi¬ 
cus as sister taxa can be seen in both ML and BA phyloge¬ 
netic inference using SAI-i. Additionally, long branches in 
the SAI-i phylograms indicate a deep divergence between 
these two species (Figs. 4 and 5). Similarly, ML phylogenet¬ 
ic inference based on protein electrophoresis data (Beer- 
li et al. 1996) revealed that P. lessonae and P. shqipericus 
shared the same ancestral node, and both were character¬ 
ized by long branches. 

Our phylogenetic inference on ND3 placed P. kurtmuel¬ 
leri in the group “ridibundus-bedriagae” sensu Plotner 
& Ohst (2001), which is the most diverse group of water 
frogs comprising at least six species (Fig. 2, Lymberakis et 
al. 2007, Akin et al. 2010a, Plotner et al. 2012). However, 
this group is supported only by a low bootstrap value in 
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1000 bp 
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Figure 3. Image of an agarose gel showing the PCR products of SAI-1 originating from three specimens each of Pelophylax kurtmuelleri 
(KK), Pelophylax shqipericus (SS) and their hybrids (KS). In P. kurtmuelleri one band of -720 bp and in P. shqipericus one of -1200 bp 
were detected. Both bands were obvious in the hybrids. 
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ML analysis (P = 70) and it is not supported by BA (Fig. 2), 
thus, adding uncertainty to the monophyly of this group. 
In line with results of the mtDNA data, the group “ridibun- 
dus-bedriagae” sensu Plotner & Ohst (2001) is not sup¬ 
ported by SAI-i (Figs. 4 and 5). 

Analyses of ND3 gene indicated low statistical support 
in ML (P=7o) and no support in BA for P. ridibundus as a 
group distinct from P. kurtmuelleri (Fig. 2). Also, branch 
lengths are short indicating a shallow genetic divergence 
between P. kurtmuelleri andP ridibundus (Fig. 2). Howev¬ 
er, sequences of P. kurtmuelleri form a supported group in 
SAI-i phylogenies (both ML and BA) (Figs. 4 and 5). 

It is also worth mentioning that the number of known 
ND3 haplotypes in P. kurtmuelleri (h = 25) is significantly 
larger than in P. ridibundus (h = 8) especially when consid¬ 
ering the larger distribution area and more analysed speci¬ 
mens of the latter. High genetic variability in species dis¬ 
tributed in southern European regions has been attribut¬ 
ed to climate oscillations during the Pleistocene (Hewitt 
2011). It can be assumed that P. kurtmuelleri had a better 
survival rate in glacial refugia across Southeast Europe 
(e.g. Lake Skadar) than P. ridibundus which probably sur¬ 
vived the Last Glacial Maximum in a yet unknown south¬ 
ern refugium from where it has recolonized Central Eu¬ 
rope. 


Hybridization between P. shqipericus 
and P. kurtmuelleri in Lake Skadar 

Analysis of SAI-i of water frogs from Lake Skadar revealed 
two fragments with a large difference in length (-720 and 
-1300 bp) (Supplementary Table Si). These fragments were 
associated to P. kurtmuelleri and P. shqipericus, respective¬ 
ly (Figs. 4 and 5). The shorter fragment was previously re¬ 
ported for P. kurtmuelleri (Hauswaldt et al. 2012) while 
the long SAI-i fragment of P. shqipericus was determined 
for the first time by Dubey & Dufresnes (2017). 

In five cases, a mixture of both types of SAI-i frag¬ 
ments was observed indicating hybridization between 
P. kurtmuelleri and P. shqipericus (Supplementary Table Si). 
In earlier studies, laboratory controlled hybridization ex¬ 
periments conducted on P. shqipericus showed its resistance 
to genome exclusion (Hotz et al. 1985). Natural hybrids of 
P. kurtmuelleri and P. shqipericus were detected, but they 
were infrequent with 13 out of 150 individuals (8.7%) and 
only two recombinant backcrosses (1.3%) (Hotz & Uzzell 
1982). These previous results are comparable to our results 
which show 5.7% hybrids among the sampled populations 
at Lake Skadar (Supplementary Table Si). 

In this research, mito-nuclear discordance, indicating 
recombinant backcrosses, was not observed in Lake Ska- 
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Figure 4. Phylogenetic tree inferred by ML analysis using SAI-1 sequences of water frogs. Newly obtained haplotypes in this study 
were marked by black rhombi. Node supports are given as bootstrap values (P), showing values > 70. 
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dar. A different situation was reported in water frogs from 
Central Europe where mito-nuclear discordance is wide¬ 
spread in hybridogenetic complexes comprising P. kl. escu- 
lentus, P. lessonae and P. ridibundus (Plotner et al. 2009). 
Hotz et al. (1985) stated that hybrids from Lake Skadar are 
not hybridogenetic and that they are only partially fertile 
(cf. Guerrini et al. 1997). This indicates that the reproduc¬ 
tion strategy and ability to induce hybridogenesis may dif¬ 
fer between - what was believed to be - P. ridibundus popu¬ 
lations from Southern Europe and those from Central and 
Western Europe (Hotz et al. 1985). As our results are in 
agreement with those of Hotz et al. (1985) and Marracci 
et al. (2011), they also support that P. kurtmuelleri is unable 
to induce hybridogenesis. The inability to reproduce hybri- 
dogenetically speaks in favour of recognising P. kurtmuel¬ 
leri as a distinct species comparing it to hybridogenetically- 
inducing but closely related P. ridibundus (cf. Holsbeek 
& Jooris 2010, Quilodran et al. 2015). However, P. kurt¬ 
muelleri potentially forms a novel hybridogenetic complex 
with P. perezi in southern France, implying that ability of 
P. kurtmuelleri may vary geographically or depends on the 
partner species involved (Dufresnes et al. (2017). 

New insights on distribution range of P. kurtmuelleri 
in Europe and conservation implications 

Quilodran et al. (2015) discussed that P. kurtmuelleri (re¬ 
ferred to as southern P. ridibundus ) represents a serious 
threat to native water frogs in Western Europe in cases of 
human-mediated translocations. This species can disrupt 
the equilibrium in hybridogenetic systems of P. kl. esculen- 
tus, which can lead to the collapse of native frogs popula¬ 
tions. The first action in the conservation of native water 


frogs should be to screen apparently non-native species in 
populations of interest. Species determination can be con¬ 
ducted relatively easily using nuclear and mtDNA genet¬ 
ic markers (Holsbeek et al. 2010, Dubey et al. 2014). For 
example, it is known that P. kurtmuelleri was introduced 
to Italy from Albania (Lanza 1962), and we found that 
the ND3 sequence “HG763874: Italy” corresponds to the 
P. kurtmuelleri lineage according to this study (haplotype 
5 in Fig. 2). In the original publication (Domeneghetti 
et al. 2013), this sequence was referred to as P. cf. ridibun¬ 
dus. Phylogenetic analyses in this study revealed that seven 
mtDNA haplotypes belonging to the P. kurtmuelleri line¬ 
age are widespread outside Southeastern Europe (Fig. 2). 
Due to a weak genetic differentiation between P. ridibun¬ 
dus and P. kurtmuelleri, these haplotypes were assigned to 
P. ridibundus or Pelophylax sp. in previous studies and not 
to P. kurtmuelleri as in this study. However, mtDNA is not 
sufficient for species affiliation of water frogs in Europe 
due to widespread mito-nuclear discordance (Plotner et 
al. 2009) and nuDNA markers such as SAI-i should be in¬ 
cluded in further analyses. 

In their native range in Albania and Montenegro, 
P. kurtmuelleri and P. shqipericus are threatened by uncon¬ 
trolled harvesting for commercial purposes (Uzzell et al. 
2009). Moreover, this is a source for illegal introduction of 
non-native species in Europe, as for instance P. kurtmuel¬ 
leri and P. shqipericus originating from Southeastern Eu¬ 
rope were reported in central Italy (Domeneghetti et al. 
2013). To cope with the increasing problem of human-me¬ 
diated translocations across Europe (Holsbeek et al. 2010, 
Dubey et al. 2014), we demonstrated the use of low-cost 
genetic determination of P. kurtmuelleri, P. shqipericus and 
their hybrids with mitochondrial and nuclear markers. The 
results of this research also provide a good basis for iden- 
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Figure 5. Phylogenetic tree inferred by BA analysis using SAI-1 sequences of water frogs. Newly obtained haplotypes in this study 
were marked by black rhombi. Node supports are given as posterior probabilities (pp), showing values > 0.9; “//” was used to present 
branch shortened to 1/8 of its original length. 
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tification of the geographic origin of P. kurtmuelleri and 
P. shqipericus genotypes in Europe in further studies. 
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